using System;
using System.Collections.Generic;
using System.Text;

namespace Barbar.MapReport.PolygonReducer {
  internal class DouglasPeucker {
    private static decimal DistanceSquared(Point u, Point v) {
      return SquaredLengthOfVector(u - v);
    }

    private static decimal SquaredLengthOfVector(Vector v) {
      return v.Dot(v);
    }

    public List<Node> SimplifyPolyline(float tolerance, List<Node> points) {
      int i, k, pv;            // misc counters
      float tol2 = tolerance * tolerance;       // tolerance squared
      Node[] buffer = new Node[points.Count];      // vertex buffer
      bool[] mk = new bool[points.Count];  // marker buffer
      Array.ForEach<bool>(mk, delegate(bool target) {
        target = false;
      });

      // STAGE 1.  Vertex Reduction within tolerance of prior vertex cluster
      buffer[0] = points[0];              // start at the beginning
      for (i = k = 1, pv = 0; i < points.Count; i++) {
        if (DistanceSquared(points[i], points[pv]) < (decimal)tol2)
          continue;
        buffer[k++] = points[i];
        pv = i;
      }
      if (pv < points.Count - 1)
        buffer[k++] = points[points.Count - 1];      // finish at the end

      // STAGE 2.  Douglas-Peucker polyline simplification
      mk[0] = true;
      mk[k - 1] = true;       // mark the first and last vertices
      DouglasPeuckerSimplificationRoutine(tolerance, buffer, 0, k - 1, mk);

      List<Node> result = new List<Node>();
      // copy marked vertices to the output simplified polyline
      for (i = 0; i < k; i++) {
        if (mk[i])
          result.Add(buffer[i]);
      }
      return result;         // m vertices in simplified polyline
    }

    /// <summary>
    /// This is the Douglas-Peucker recursive simplification routine
    /// It just marks vertices that are part of the simplified polyline
    /// for approximating the polyline subchain v[j] to v[k].
    /// </summary>
    /// <param name="tolerance">approximation tolerance</param>
    /// <param name="v">polyline array of vertex points</param>
    /// <param name="j">indices for the subchain from v[j]</param>
    /// <param name="k">indices for the subchain to v[k]</param>
    /// <param name="mk">mk[] = array of markers matching vertex array v[]</param>
    void DouglasPeuckerSimplificationRoutine(float tolerance, Point[] v, int j, int k, bool[] mk) {
      if (k <= j + 1) // there is nothing to simplify
        return;

      // check for adequate approximation by segment S from v[j] to v[k]
      int maxi = j;          // index of vertex farthest from S
      decimal maxd2 = 0;         // distance squared of farthest vertex
      decimal tol2 = (decimal)(tolerance * tolerance);  // tolerance squared
      Segment S = new Segment(v[j], v[k]);  // segment from v[j] to v[k]
      Vector u = S.P1 - S.P0;   // segment direction vector
      decimal cu = u.Dot(u);     // segment length squared

      // test each vertex v[i] for max distance from S
      // compute using the Feb 2001 Algorithm's dist_Point_to_Segment()
      // Note: this works in any dimension (2D, 3D, ...)
      Vector w;
      Point Pb;                // base of perpendicular from v[i] to S
      decimal b, cw, dv2;        // dv2 = distance v[i] to S squared

      for (int i = j + 1; i < k; i++) {
        // compute distance squared
        w = v[i] - S.P0;
        cw = w.Dot(u);
        if (cw <= 0)
          dv2 = DistanceSquared(v[i], S.P0);
        else if (cu <= cw)
          dv2 = DistanceSquared(v[i], S.P1);
        else {
          b = cw / cu;
          Pb = S.P0 + u * b;
          dv2 = DistanceSquared(v[i], Pb);
        }
        // test with current max distance squared
        if (dv2 <= maxd2)
          continue;
        // v[i] is a new max vertex
        maxi = i;
        maxd2 = dv2;
      }
      if (maxd2 > tol2)        // error is worse than the tolerance
    {
        // split the polyline at the farthest vertex from S
        mk[maxi] = true;      // mark v[maxi] for the simplified polyline
        // recursively simplify the two subpolylines at v[maxi]
        DouglasPeuckerSimplificationRoutine(tolerance, v, j, maxi, mk);  // polyline v[j] to v[maxi]
        DouglasPeuckerSimplificationRoutine(tolerance, v, maxi, k, mk);  // polyline v[maxi] to v[k]
      }
      // else the approximation is OK, so ignore intermediate vertices
      return;
    }
  }
}
